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We study the propagation of a massive vector or Proca field on the Schwarzschild spacetime. The 
field equations are reduced to a one-dimensional wave equation for the odd-parity part of the field 
and two coupled equations for the even-parity part of the field. We use numerical techniques based 
on solving (scalar or matrix- valued) three-term recurrence relations to compute the spectra of both 
quasinormal modes and quasi-bound states, which have no massless analogue, complemented in the 
latter case by a forward-integration method. We study the radial equations analytically in both the 
near-horizon and far-field regions and use a matching procedure to compute the associated spectra 
in the small-mass limit. Finally, we comment on extending our results to the Kerr geometry and its 
phenomenological relevance for hidden photons arising e.g. in string theory compactifications. 
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I. INTRODUCTION 

The propagation of massless fields and linearized met- 
ric perturbations on black hole spacetimes has been ex- 
tensively studied over the last half-century, following the 
foundational work of Regge and Wheeler [l|. Key mo- 
tivations for such studies range from questions of black 
hole stability @ , no-hair theorems 01 , Hawking radiation 
[Hi, scattering and absorption 0-0]; an d gravitational- 
wave dynamics In the 1970s, it was established that 
the response of a black hole to a generic initial pertur- 
bation is characterized by a spectrum of damped reso- 
nances Igl, c alled quasinormal modes (for recent reviews 
see e.g. (10L [Tlj ) and that, after a period of quasinor- 
mal ringing, the perturbing field subsequently undergoes 
power- law decay [l2[ . 

The propagation of massive fields, such as the electron, 
neutrino and baryonic fields, on black hole spacetimes 
has also received some attention [l3l - l32j . If the field has 
mass fi and the black hole has mass M, the key dimen- 
sionless quantity is the mass coupling Mfi = G Mfi,/ He = 
Mfi/rrip where mp is the Planck mass (henceforth we set 
G = c = tl = 1 ) . For astrophysical black holes and known 
elementary fields we expect Mfi to be large. For exam- 
ple, for a (hypothetical) neutrino mass fi u ~ 0.02 cV, one 
obtains Mfi = 1.5 x 10 8 for solar-mass black holes. We 
therefore expect massive fields to behave very differently 
from massless fields in most astrophysical scenarios, at 
least for multipoles I < Mfi. 

Time domain studies p7l - l29l . [32[ have also revealed in- 
teresting differences in the propagation of massive and 
massless fields. First, quasinormal ringing is generally 
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less easily excited, and the frequencies and decay times 
depend on Mfi. Second, power-law decay is replaced 
with a slowly-decaying, oscillating phase (reminiscent of 
Zitterbcwegung) in which the frequency of oscillation is 
given by the mass of the field, rather than the mass of 
the black hole. It has been shown that, at late times 
fit ^s> l/(Mfj,) 2 , all multipoles exhibit a 'universal' oscil- 
latory decay, with a spin- independent envelope ~ i~ 5 / 6 
[2714291 . HH and frequency fi. Third, massive fields with 
spin exhibit additional longitudinal radiative degrees of 
freedom, and low multipoles I < s of massive fields (for 
example the vector field monopole) have a radiative char- 
acter, rather than the static character apparent in the 
massless case. It has long been known that a black hole 
cannot support massive 'hair', and consequently does not 
have a well-defined baryon number Q . 

From the point of view of phenomcnonology, the most 
interesting scenarios arise when Mfi ~ 1, a condition 
that may be attained either for light primordial black 
holes postulated to be produced during inflation [33l - [35j 
or 'ultralight' exotic particles found in beyond the Stan- 
dard Model scenarios. These scenarios are particularly 
significant in the context of the well-established super- 
radiant instability of rotating black holes [IH 
H3| subject to massive bosonic bound states. For the 
massive scalar field, numerical studies have found (42r - 
l44l . [46| that the maximum growth rate of the instabil- 
ity is approximately r _1 ~ 1.5 x 10~ 7 (GA//c 3 ) -1 for 
a rapidly-rotating black hole at a ~ 0.99M, and it oc- 
curs for a coupling Mfi ~ 0.42. Hence, this process is 
relevant in astrophysical environments for particles with 
masses fi < 10~ 10 eV. Ultra-light particles have been pro- 
posed to be ubiquitous in string theory compactifications 
[48j |. where the several axionlike fields acquire masses 
exclusively through nonperturbative string instanton ef- 
fects, making them exponentially sensitive to the size of 
the associated compact cycles. This leads to a generic 
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landscape of ultralight axions, known as the 'string axi- 
verse' |48l - [5G |. populating all mass scales possibly down 
to the present Hubble scale, Ho ~ 10 -33 eV. Recently, 
there have been suggestions that ultralight bosonic fields 
(e.g. /i ~ 10 -22 eV) may play the role of dark matter 
in galactic halos [51(. It is possible that rotating black 
holes may provide a unique probe of the ultralow energy 
spectrum of string theory compactifications, with an ex- 
tremely rich phenomenology associated with the forma- 
tion of bosonic superradiant bound states around astro- 
physical black holes, such as the emission of gravitational 
waves, gaps in the mass-spin black hole Regge spectrum 
and 'bosenova-like' particle bursts [49j . 

Massive hidden U(l) vector fields are also a generic 
feature of Bcyond-Standard-Modcl scenarios and, in par- 
ticular, string theory compactifications, arising in the lat- 
ter case from a variety of sources, such as broken non- 
Abelian orbifolds in heterotic compactifications, and D- 
brane configurations and bulk Ramond-Ramond fields 
in type II string theories p52l - [54l ]. The nature of U(l) 
vector masses is, however, inherently different from the 
axion case, where the associated shift symmetry (remi- 
niscent of an underlying U(l) symmetry) is only broken 
nonperturbatively. U(l) gauge symmetries may, on the 
other hand, be broken by perturbative effects following 
the Higgs mechanism or by Stueckelberg-type couplings. 
In the former case, the resulting vector masses are gener- 
ically determined by the soft-supersymmetry breaking 
terms in the hidden sector and hence are quite sensi- 
tive to the mechanism mediating supersymmctry break- 
ing to this particular sector, which may yield a broad 
range of possible values. Light U(l) particles can nev- 
ertheless be obtained in this case for hyperweak hidden 
sector gauge couplings, which naturally arise e.g. from 
volume suppression if the hidden sector lives on D-branes 
wrapping a large compact cycle. Alternatively, one may 
envisage scenarios where the string coupling is actually 
extremely small in the overall compactification (within 
the limits of the 'gravity as the weakest force' conjec- 
ture [Hj]), whereas the visible (Standard Model) sector 
wraps a collapsed cycle that makes the associated gauge 
couplings unnaturally large f60j | . 

Given such motivations, it is natural to speculate on 
the possibility that black hole physics, and in particular, 
the mechanism of unstable bosonic bound states driven 
by superradiance, can be used to probe the existence of 
ultralight hidden vector particles. For example, super- 
massive black holes, with masses ~ 10 5 — 10 10 M so i ar are 
believed to be found in the centre of most active galax- 
ies. Unstable 'bosonic clouds' may form in the vicin- 
ity of rotating supermassive black holes, if particles with 
masses in the range ~ 10~ 16 — 10~ 21 cV exist. (Note 
that a reasonable upper bound on the mass of the pho- 
ton is m 7 < 1CP 18 eV (EH). Stellar-mass black holes, 
formed in supernovae, are likely to be even more abun- 
dant; such black holes are sensitive to particles with 
masses ~ 10~ 10 — 1CP 12 eV. Since the amplification mech- 
anism relies on purely gravitational physics, it is largely 



independent of the kinetic mixing between the hidden 
and visible photons which is the basis for other types of 
analysis (see e.g. [H]). 

There are also avenues to explore within the context 
of Standard Model physics. For example, in astrophys- 
ical environments, the electromagnetic field can acquire 
a small 'effective mass' via the Anderson-Higgs mecha- 
nism, induced by the accreting plasma surrounding the 
black hole [13, HH ■ This raises the possibility that super- 
radiant instabilities may even play a role during black 
hole accretion. As a step in this direction, the effect of 
strong magnetic fields on the instability timescale of the 
massive scalar field was considered in [451 ] . 

The study of massive vector field perturbations in ro- 
tating black hole spacetimes poses a challenging problem 
from both the analytical and numerical points of view, 
due to the apparent nonseparability of the Proca equa- 
tions on the Kerr spacetime. Nonseparability is particu- 
larly frustrating because all other massless and massive 
fields (at least scalar and Dirac) admit separable solu- 
tions. Since it seems that any exploration of the phe- 
nomenology of the Proca field on the Kerr spacetime will 
require a careful analysis of partial (rather than ordi- 
nary) differential equations, it is essential that we first 
establish a solid understanding of the Proca field on the 
Schwarzschild spacetime, where the equations are sepa- 
rable, which is the purpose of this work. Of course, since 
the Schwarzschild spacetime does not exhibit superradi- 
ance, we do not expect to find any instability in this case. 

Massive vector field perturbations on spherically- 
symmetric spacetimes have received some attention in 
the literature [13, HH IH, [59|, although comparatively 
less than other massive and massless perturbations. The 
problem of solving the massive vector field or Proca 
equations on the Schwarzschild spacetime is in itself 
quite challenging, given the coupling between the differ- 
ent components of the field that results from the bro- 
ken gauge invariance and the additional longitudinal de- 
gree of freedom. Previous studies have examined the 
monopole mode [U l32ll . described by a single equation, 
and large multipoles [T7j, which can be studied using 
WKB methods. Very recently, the Proca field in higher- 
dimensional spherically symmetric spacetimes was exam- 
ined in the context of Hawking radiation [59j . Here, we 
analyze Proca field perturbations on the Schwarzschild 
spacetime for generic multipoles, using both numeri- 
cal and analytical techniques to focus on the spectrum 
of quasinormal (QN) modes and (quasi-)bound states. 
These spectra reveal several peculiar features that are 
related to both the massive and the higher-spin nature 
of the perturbations. To interpret our findings, we then 
compare our results with those for both the massive 
scalar and the (massless) electromagnetic fields. 

This work is organized as follows. In Sec. |H] we out- 
line the properties of the Schwarzschild spacetime and 
formulate the Proca equations describing massive vec- 
tor field perturbations on this geometry In Sec. IIIBl we 
separate the radial and angular parts of the spin-1 wave- 
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function using vector spherical harmonics, and we show 
that odd- and even-parity parts are completely decou- 
pled. In Sec. Ill CI we explore the relationship between 
the massless limit of the Proca field and the electromag- 
netic field. In Sec. IIIII we define the massive spin-1 QN 
modes and bound states, and describe numerical methods 
for computing their spectra, based on continued-fraction 
and forward-integration techniques. We present a selec- 
tion of numerical results for the spectra in Sec. [IVJ In 
Sec. El we develop analytical techniques to study the ra- 
dial wavefunctions in the near- and far-field regions, and 
apply a functional matching procedure to obtain the as- 
sociated spectra in the limit Mfj, -C 1. We conclude with 
a summary of our main results in Sec. IVI| where we also 
discuss possible extensions and future work on this topic. 



II. ANALYSIS 
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where Y lm = Y lm (9,4>) denote the ordinary (scalar) 
spherical harmonics and sg = s'm9. Note that the four 
vector spherical harmonics are defined in an analogous 
way to the ten tensor spherical harmonics used for de- 
scribing gravitational perturbations in the Lorenz gauge 
[60| . These harmonics satisfy the orthogonality condi- 
tions 



(8) 



A. Spacetime and field equations 

The Schwarzchild spacetime is described by coordi 
nates x M = {t, r, 9, <fi} and the line element 

j-2 , i-ljj , 2 ( jq2 i • 2 a j ±2\ 



where we defined r]^ v = diag[l, f 2 , 1/r 2 , l/(r 2 sin 2 9)] and 
dVL = sm9d9d(f>. 

Let us briefly consider the transformation properties 
of these harmonics under parity inversion, x — > — x, i.e. 
9 — > it — 9 and <j) ~ * § + under which the spatial 
components of a vector field transform as [A r , Ag, A$] —> 
ds 2 = ~f(r)dt 2 + r l dr 2 + r 2 (d9 2 + sin 2 9d(j) 2 ) , (1) [_A r , +A g , -A+\. The first three harmonics (i = 1, 2, 3) 

in Eq. (j4|) then pick up a sign of (—1)' under parity in- 
version, whereas the last (i = 4) harmonic picks up the 
opposite sign (— 1)' +1 . It is standard to call the former 
the 'even-parity' or 'electric' modes and the latter the 
'odd-parity' or 'magnetic' modes. 

We may now decompose the vector potential in this 
basis: 



where f(r) = 1 — 2M/r. The exterior region (?* > th = 
2M) is stationary and static; in other words there exists 
a timclikc Killing vector (^ = Sq ) which is hypersurfacc- 
orthogonal. The spacetime is Ricci-flat, i.e. R = and 
R[iv = R^ 0. 

The Proca field strength is defined in terms of the 



(physical) vector potential A v by F^ v = A v;IJi — A^ v , ^ 4 

where ift denotes a covariant derivative with respect to A^(t,r,9,(j)) = - Cj u l $(t, r)Z^ lm (6, 4>) , (9) 



The field equations in vacuum are simply 



F» v . ll = l x 2 A v 



(2) 



It follows that F^ v 



and hence the Lorenz condition 



A?„ = is simply a consequence of the field equations. In 
other words, there is no gauge freedom in the Proca field, 
which thus describes three physical degrees of freedom. 
Hence, the field equations can be written as 



where c\ = c 2 = 1, C3 = C4 = [1(1 + l)]~ x / 2 . This ansatz 
is sufficient to separate the vector field equations ^ into 
a set of four second-order partial differential equations in 
r and t, 



A"' 11 ,, - [i 2 A v = 



(3) 



with the supplementary condition A^ — 0. In the fol- 
lowing sections, we will often take M = 1. 



B. Separation of variables 

To separate the angular part of the vector potential 
components, we introduce a basis of four vector spherical 
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Here u = v! = -pr and the tortoise coordinate r + is 

at 5 or* a * 

defined via dr* = f~ 1 dr. The differential operator 2? 2 is 
given by 



'dt 2 + drl * 



1(1 + 1) 



(14) 
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The fourth equation (jl3|) . which describes the odd-parity 
sector, is completely decoupled from the first three equa- 
tions (|T0|) -([T2 j) . which describe the even-parity sector. 
The even-parity equations must be supplemented by the 
Lorenz condition, 



Hi) 



*(2) 



'(2) 



'(3) 



)=0. 



(15) 



This condition may be used to reduce the even-parity 
system to a pair of coupled differential equations. We 
may then replace Eq. ([TT]) with 



(2) 



2/ 



1 - 



l W (2) 



'(3) 



(16) 



and note that Eqs. (| 12f) and (fTBj) form a closed system. 
Hence, we have obtained one decoupled wave equation 
for the odd-parity part of the vector potential and two 
coupled wave equations for the even-parity modes. 



C. Proca field and the electromagnetic limit 

Let us consider the relation between the Maxwell field 
and the massless limit of the Proca field. 



1. The monopole mode 

In the special case of the monopole mode (I = 0, even- 
parity), for which only the first two (i = 1 and 2) har- 
monics are defined, we obtain a single decoupled equation 

'2(r-3) 



dt 2 dr 2 



'(2) 







(17) 



This is precisely the monopole equation investigated by 
Konoplya (see Eq. (15) in [3l[). The Lorenz condition 
([T5|) implies that = £<9 r (rU(2)). 

In the massless limit (u — > 0) we may apply a 'Chan- 
drasekhar transform' [6l[ [62j to the monopole equation 
to show that iim satisfies a scalar wave equation, 

dt 2 dri r 3 







(18) 



In the electromagnetic limit, this degree of freedom is not 
physical and can be removed by a gauge transformation. 
In other words, in clcctromagnctism the monopole part 
of the field does not have radiative degrees of freedom 
and, in the Lorenz gauge, can be written as the sum 
of a gauge mode and a static field, i.e. 



A„ 



'-6° 



(19) 



where \ is arL arbitrary scalar field satisfying = 0. 

In the Proca case (fj, ^ 0), the static mode in Eq. (fT9")) 
is not a solution of the equation, whereas the radiative 
degree of freedom is physical. In the small-mass limit, 
we then expect the spectrum of the monopole part of the 
Proca field, governed by (|17[) . to approach the spectrum 
of the I — mode of a scalar field. 



2. Massless limit of odd-parity modes: Price 's equation 



In the massless limit, the odd-parity equation (|13p re- 
duces to Price's equation fl2j . 



V 2 {r 2 4>i) = 



(20) 



which governs the evolution of the Maxwell scalar of spin- 
weight zero <f>\. In fact, inserting the ansatz in Eq. (jH]) for 
A p into the definition of <f>\ = F iiv (l^n 1 ' + m^m"), where 
n M , m M and m M are the null vectors of the Kinnersley 
tetrad 1631. leads to 



1(1 + 1) 



(21) 



3. Massless limit of even-parity I > modes 



The pair of coupled wave equations for the even-parity 
sector, Eqs. (JT3J) and (jTHJ) , can be combined into a single 
fourth-order equation. In the massless case, fx = 0, it 
is notable that the fourth-order equation for it( 3 ) can be 
factorized as follows: 



j_,j=J,S=l 

r p U RW 



1^1^=0 
RW 



l (3)) 







(22) 



where we have defined the generalized 'Regge- Wheeler' 
operator T> l ^ w as 

d 2 d 2 



V RW~ dt 2 + dr 2 f 



1(1 + 1) 2(1 -s 2 ) 



(23) 



The equation for it( 2 ) can also be factorized, although the 
operators do not have such a simple form. 



It can be shown that if : 



*(3) 



satisfies the 'inner' scalar 



wave equation in Eq. (1221) . i.e. T)^ w (ru^) = 0, then the 
corresponding vector potential is pure-gauge, i.e. 



A„ 



X, 



where 



X 



i(3) 



1(1 + 1) 



y Zm (M) ■ (24) 



In other words, in the electromagnetic case, the 'inner' 
degree of freedom has no physical significance since it 
can be completely removed by a gauge transformation 
A^ —> A^ — x,n- However, for /i ^ there is no gauge 
freedom and so this degree of freedom becomes physical. 
We expect its spectrum to approach that of a scalar wave 
in the low-mass limit, as for the monopole mode. 

The 'outer' wave equation in Eq. (|2"2"|) for the quantity 



ip = / _1 2?^,'^ (ru(3)) has physical significance even in 
the massless case. It is straightforward to show that ip 
can be written as 



^ = "(3) ~ U (2) 



(25) 



Working in reverse, inserting Eq. (j2"5j) into Eqs. (|12l) and 
([Ttll) leads to a 'vector' (s = 1) wave equation for ip, 



= V 2 t/J = 



(26) 
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We can conclude that, in the electromagnetic limit, 
the physically-meaningful even-parity and odd-parity de- 
grees of freedom are described by the same dynamical 
equation [see Eqs. (fTB")) and (|2"6")l ] and hence have the same 
spectrum. However, this degeneracy is broken if the field 
has mass; in addition, the gauge degrees of freedom (in 
the monopole and in higher modes) acquire physical sig- 
nificance in the Proca case. In the following sections, we 
investigate the richer structure of the Proca spectrum. 

III. FREQUENCY SPECTRA: QUASINORMAL 
MODES AND BOUND STATES 

In this section, we investigate the spectrum of charac- 
teristic modes of the Proca equation on the Schwarzschild 
spacetime. We start with a frequency-domain represen- 
tation, where 

«^(t,r) = «^( W ,r)e- &rt (27) 

and oj may take complex values. Note that Im(w) < 
corresponds to exponential decay, while Im(w) > leads 
to exponential growth, such that we expect to find only 
the former behaviour in the Schwarzschild spacetime. 

Physical modes on a black hole spacetime must be 
purely ingoing at the horizon from the point of view of a 
local observer. The ingoing condition corresponds to the 
horizon boundary condition 

« (i) (w J r)~e-*- r - , (28) 

as r* —> — oo. In the asymptotically flat region, the solu- 
tion resembles 

u (i) (w, r) ~ B {i) (u)e- kr * + (u,)e +kr * , (29) 

as r* — > +oo, where B^(ui), C^(oj) are complex coeffi- 
cients and we define k = ji 2 — up such that Rc(fc) > 0. 
There are two types of special mode that we may con- 
sider: (i) the QN modes, defined by B^(oj) = 0, which 
asymptotically resemble purely outgoing waves; (ii) the 
quasibound states, defined by C(j)(w) = 0, which are 
spatially localized within the vicinity of the black hole, 
i.e. decay exponentially away from the black hole. In 
either case, imposing an asymptotic boundary condi- 
tion generates a discrete spectrum of allowed frequen- 
cies, {w;^(M/x)}. Here, each QN mode or bound-state 
frequency is labeled by its angular momentum number 
I > 0, overtone number n > 0, and 'polarization' num- 
bers A, which we further describe below. 



recurrence relations for series coefficients, which may be 
found by solving a continued- fraction equation (64|. An 
appropriate ansatz is 

(oj, r) = f-^ r - v e" r £ a« [f(r)} n , (30) 

n 

where v = (oj 2 — q 2 )/q. To seek QN frequencies, we set 
q = k [i.e. q = \J /i 2 — oj 2 and thus Rc(q) > 0]. On the 
other hand, for bound-state frequencies we set q = —k 
[i.e. Re(q) < 0]. Inserting Eq. (|30| into the governing 
equations ffl2|) . (|T3|) . (fT6|) and (fTTj) then leads to three- 
term relations of the form 

a ai + /3 a = (31) 
a n a n+1 + f3 n a n + 7„a„_i =0, n > . (32) 

Where the equations are decoupled, i.e. for the odd-parity 
and monopole cases, a„, /3 n , j n are scalar quantities, and 
the QN or bound-state frequencies are those for which 

a Qn-l7n _ 

Pn r, a n -27n-i ~ 

P ™" 1 0„^ 2 -a„_37„-2/... 

ID Qn + l7n + 2 ' {^J) 

0n+2— a„+27n+3/--- 

Methods for solving continued-fractions have been de- 
scribed in detail elsewhere [(34[ and this method is exten- 
sively used (see e.g. Sec. IIIH in [ll| or J65|), so that we 
will give no further details here. Where the equations 
are not decoupled, i.e. for even-parity modes, a n etc. are 
matrix-valued. This case is discussed in Sec. MI A 31 



1. Monopole mode 

The monopole mode I = is governed by E q. (117|) . 
The QN mode spectrum was found by Konoplya [3I( via 
the continued-fraction method. We obtain the following 
coefficients: 

a n = (n + l)(n+l-4iu>) , (34) 
P n = -2n 2 - 2q~ 1 [Aiojq + 3q 2 - lo 2 + q) n 
+q~ 1 (4iw 3 + 12qui 2 - Yliuiq 2 - 4q 3 + 
+q-uj 2 + Aiuq + 3q 2 ) , (35) 
7„ = n 2 - q^ 1 (4:iujq + 2q 2 - 2u 2 ) n + 

+q-' 2 (u - iqf - 4 . (36) 



A. Continued-fraction method 

As described in [44[ , either frequency spectrum may be 
found using a method which relies on the fact that, with 
a suitable ansatz, the desired solutions of the differential 
equations correspond to minimal solutions of three-term 



2. Odd-parity modes 

The odd-parity modes satisfy the decoupled wave 
equation (fT3"|) which reduces to the s = 1 Regge- Wheeler 
(RW) equation T^^^u^ = in the massless limit. In 
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2 2 \ 

uj — q)n 



this case, we obtain a n as above and 
/3 n = -2n 2 + 2q~ 1 (Aituq + 'Sq 2 

+q- 1 [3q 2 + {12q - 1)oj 2 - 4q 3 - ql(l + 1) - 
-\2iujq 2 + Aiujq + Aiuj 3 ], (37) 
7„ = n 2 — 2q~ 1 (q 2 — uj 2 + 2iojq) n — 

- q- 2 [- iiuq 3 + 4:iuj 3 q+6uj 2 q 2 -q i -uj 4 + q 2 ] . (38) 

3. Even-parity modes 

The even-parity modes satisfy a pair of coupled differ- 
ential equations, Eqs. (|12j) and (|16|). Inserting the ansatz 
(|30| into these equations then leads to a matrix-valued 
three-term recurrence relation, 



anU n+ i 
with a vectorial coefficient U„ 



a Ui + A,U = , (39) 
/3„U„ + 7nU„_! =0, n > , (40) 

/„( 2 )\ 



,(3) 



and matrices: 



In 



av, 

7« - 3 3 
7. 



/3„ = 



A, + 1 -1 
21(1 + 1) /3 n 



(41) 



where 



(n + l)(n + 1 - Aioj) , 

~6<7 2 + 8igcj-2w 2 



-2n 2 



1 



1(1 + !) + -[- 12iujq 2 + 3q 2 - 4q 3 



- uj 2 + \2quj 2 + Aiujq + Aiui 2 
(uj — iq) 2 



(42) 



In 

The matrix- valued three-term recurrence relation can 
be solved using matrix- valued continued fractions [66| . 
We search for roots of the equation MU = 0, where 

M = fa - a Q [f3i - on (/3 2 + a 2 R^) 72] ^ H , (43) 
with U n+ i = R+U„ and 

R+ = - (/3 n+ i + a„+iR+ +1 ) ^ 7„+i . (44) 

For nontrivial solutions Uo, one must hence solve 

det|M|=0. (45) 

In nondegenerate cases, we expect to find two indepen- 
dent solutions, with distinct eigenvectors and Ug . 
To distinguish between the solutions, we define the quan- 
tity 



lim 



E«i 3) /E4 2) 

n=0 ii=0 



which we will loosely refer to as the polarization. 



(46) 



B. Forward-integration method 

An alternative numerical method can be used to com- 
pute the spectrum of bound states, based on the expected 
convergence of the solutions at infinity. In the vicinity 
of the horizon, the expansion in Eq. (j30j) for the radial 
functions satisfying ingoing boundary conditions may be 
written as 



i(*)( w > r ) = ( r - r H) 



OO 

n=0 



■)n(l 



(47) 



The series coefficients for n > 1 can be easily found, 
as a function of the leading coefficient 6(i)o, by substitut- 
ing this ansatz into the radial equations (|12p . (|13p and 
(|T6|) and using a symbolic algebra package, e.g. Mathc- 
matica. We may then use Eq. (jT7|) as an initial condition 
close to the horizon to numerically integrate the radial 
equations up to the far-field region. 

For the odd-parity modes, we expect this to lead to 
an asymptotic solution of the form Eq. (|29p . i.e. a linear 
combination of an exponentially divergent and an expo- 
nentially convergent term, with the former vanishing for 
bound states. We may thus determine the bound-state 
spectrum by setting 6(4)0 = 1 an d minimizing the result- 
ing U(4)(w,r) for an arbitrarily large distance r S> rn in 
the complex w-plane. 

This method can be extended for the even-parity 
modes, where the mode equations are coupled. In this 
case, we may obtain a family of numerical solutions at in- 
finity parametrized by the unknown leading coefficients 
(6(2)o> ^(3)0)) with bound states corresponding to partic- 
ular values of the latter that lead to pure exponentially 
decaying solutions. The associated spectrum may then 
be computed by choosing a suitable basis for the space of 
initial coefficients, for example (6(2)0^(3)0) = (1,0) and 
(0, 1), and defining a 2 x 2 matrix of solutions 



S(w,r) 



t (1,0)/ \ (0,1)/ x 
u (2) K r ) U ( 2 ) K r ) 



1,0)/ \ (0,1)/ \ 
3) ( w ' r ) "(3) ( W ' r ) 



(48) 



"(3) 



The particular linear combinations of the near-horizon 
solutions corresponding to bound states will thus corre- 
spond to the kernel of S evaluated at r — > 00. In practice, 
this corresponds to minimizing det |S| in the complex to- 
plane at an arbitrarily large distance. Notice that, for 
each minimum, only the eigenvector associated with an 
asymptotically vanishing eigenvalue will correspond to 
a physical solution, whereas the remaining eigenstate is 
unphysical and yields an exponentially large eigenvalue. 
This allows one to reconstruct the radial wave functions 
in each case and determine the associated polarization, 
as defined in Eq. (|4*p| . 

Although this method cannot be applied to compute 
the QN mode spectrum, where the purely divergent na- 
ture of the solutions is hard to determine numerically, it 
can be more easily implemented e.g. with Mathematica, 
where the lengthy algebraic expressions resulting from 
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the multiple matrix inversions required by the continued- 
fraction method are rather difficult to minimize. Further- 
more, forward-integration provides an independent check 
of the numerical results obtained with the latter method, 
making our analysis more robust. 



IV. NUMERICAL RESULTS 
A. Quasinormal modes 

Fig. Q] shows the effect of mass on the quasinormal fre- 
quencies of the \ow-l, n modes of the Proca field. For a 
given angular momentum number I and overtone number 
n, there are two even-parity modes and one odd-parity 
mode. The even-parity modes may be distinguished by 
their behaviour in the massless limit, as discussed in 
Sec. Ill CI In this limit, the spectrum of 'scalar' modes 
(which are unphysical pure-gauge modes in electromag- 
netism) reduces to the spectrum of a scalar field. In 
the same limit, the 'vector' even-parity and odd-parity 
modes are degenerate, with the frequencies of the elec- 
tromagnetic field. The field mass breaks this degeneracy. 

Quasinormal mode spectrum of Proca field on Schwarzschild spacetime 
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FIG. 1: Quasinormal mode frequencies of the Proca field, for 
I — (monopole), / = 1 (dipole) modes, I — 2 (quadrupole) 
modes, for a range of field masses Mfj, = 0, 0.04, . . . , 0.2. The 
fundamental (n = 0) and first overtones (n = 1) are shown. 
For a given I, n there are two even-parity modes, and one 
odd-parity mode. In the massless limit, the 'scalar' even- 
parity mode has the same QN frequency as the scalar (s = 0) 
field, whereas the 'vector' even-parity and odd-parity modes 
have the same QN frequency as the electromagnetic field. 

The plot in Fig. [T] illustrates the effect of the field mass 
upon the QN spectrum of low-Z, low-n modes. As ex- 
pected, this is greatest for the lowest modes. As previ- 
ously observed by Konoplya [3l|, the decay rate of the 
monopole mode decreases substantially as the mass cou- 
pling M n is increased. However, the physical relevance 
of the QN mode also diminishes as the mass increases, as 
the height of the effective potential barrier decreases. 



In the Proca case (fi ^ 0), the two even-parity QN 
modes ('scalar' and 'vector') have distinct frequencies 
and polarization states. In Fig. [3J we examine the po- 
larization of the modes at large distances, by plotting V 
as defined in Eq. (|46[) function of M/x, for i = 1,2,3 
and for the fundamental mode (n = 0) and the first over- 
tone (n = 1). In the massless limit, V approaches 
for 'scalar' modes and 1 for 'vector' modes. Addition- 
ally, one observes that V varies smoothly with Mfi. In 
theory, a measurement of V and wz„ would allow one to 
independently deduce the mass of the black hole M and 
the mass of the field \i. In practice, however, detecting 
QN ringing from a Proca field is unlikely to be possible 
in the foreseeable future. 

Polarization for Quasinormal Modes of Even-Parity Fieid for Mass Couplings Mu = .. 0.2 
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FIG. 2: Polarization state of even-parity quasinormal modes. 
The plot shows the complex number V, i.e. the ratio u^)/u(2) 
far from the black hole [defined in Eq. (|46|) ]. as a function of 
the mass coupling Mfj, = . . . 0.2, for I = 1, 2, 3 and n = 0, 1. 
For scalar (vector) modes, V — > (V — ► I) as Mfi — ¥ 0. 



B. Bound states 



As discussed in [13|, Ll4, ll6|, Ha l20|, l21|, [2^, l26|, [40|, , 

a massive field may be localized in the vicinity of a 
black hole in (quasi-)bound states with complex frequen- 
cies. In this section we present a selection of numeri- 
cal results for the Proca-field bound-state spectrum. We 
have verified that the results obtained via the continued- 
fraction method (Sec. IIII A[) are in excellent agreement 
with those obtained with the forward-integration tech- 
nique (Sec. IIII Bp . 

The bound states of the Proca field were previously 
considered in Ref. [TtJ for large multipoles, where it was 
shown that in the limit M/i — > 0, the spectrum is hydro- 
genic, i.e. [70[ 



Rcfw/V) w 1 



(M/i) : 



(49) 



2N 2 ' 

where N — j + 1 + n and j — I + S is the total angular 
momentum of the state as measured by an asymptotic 
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observer, with spin projection S = 0,±1. Our results 
are fully consistent with this if S = +1 for the monopole 
mode, in agreement with the rules for addition of angular 
momenta, such that \l — 1| < j < I + 1 for spin-1 fields. 

In Fig. |3] we show the bound-state frequency spectrum 
ui/fj> as a function of the mass coupling M fx, for the lowest 
modes I = 0, 1. For a given I and n, there are three types 
of mode: (i) odd parity, S — 0, (ii) even parity, S = +1 
and (iii) even parity, S — — 1, the monopole being a type 
(ii) mode. As predicted by Eq. ((49]), the lowest-energy 
mode is I = 1, S = — 1. 

Energy levels of bound slates of Praca field on Schwarzschild spacetime 
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Energy levels of bound states of Proca field on Schwarzschild spacetime 
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FIG. 3: Bound state levels of the Proca field on the 
Schwarzschild spacetime. The upper plot shows the real part 
of the frequency Re(w//u) as a function of the mass coupling 
M/j, and the lower plot shows (the negative of) the imaginary 
part Im(cj/^) on a logarithmic scale. The modes are labeled 
by their angular momentum number I, overtone number n, 
spin projection S and parity (odd or even). 

The imaginary part of the frequency, which sets the de- 
cay rate of the mode, increases monotonically with M/x. 
The lower plot of Fig. [3] shows that, in the small- Mfx 
regime, there is a power-law dependence, Im(u;//x) oc 
— (Mll) v , where 77 depends on I and spin projection S, 
and the constant of proportionality depends on the over- 
tone number. From the numerical data we infer that 



77 = 4/ + 25 + 5. 



(50) 



For example, Fig.@]shows the exponent 77 estimated from 
the numerical data for Mfi <C 1, clearly showing that the 



modes I = L, S = +1 and I = L + 1, S = — 1 have the 
same exponent rj. Our data is not in agreement with the 
results found in [17j . suggesting that the latter analysis 
is not applicable to the lowest multipoles. 



Power-law index ri in decay rate lm((o/u) oc -(Mfi) 11 : Numerical data 
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FIG. 4: Numerical data for the exponent r\ in the power- 
law relationship Im(w//i) oc —(Mjj,) v which determines the 
decay rate of the quasi-bound states in the small-coupling 
regime Mfi <C I- Here I — 0, 1, . . . is the angular momentum 
number, and S G { — 1,0, +1} is the spin projection of the 
state in the large-r regime. The data shown were obtained by 
numerically evaluating the function 77 = — ^-f^pj^p^- Note 
that numerical evaluation becomes increasingly difficult in the 
small-M/j regime, where \lm(to / fj,)\ is tiny (< 10 -12 ). The 
data strongly suggests that modes I = L, S = +1 and I — L + 
1,S — —1 share the same exponent 77. The data is consistent 
with Eq. (|50|) . which implies that, in the limit Mp — > 0, the 
exponent tends to 77 = 7 (I — 0, S = +1 and 1 = 1, S = —1), 
77 = 9 (I = 1, S = 0) and 77 = 11 (I = 1, S = +1 and 
1 = 2, S= -1). 



The pair of even-parity modes may be distinguished 
by examining the polarization V defined in Eq. (|46p . For 
S = +1 modes we find P — > — I as Mfj, — > 0, and for 
S = — 1 modes we find V — > I + 1 in the same limit, as 
shown in Fig. [S] 

Fig. [HI shows the comparison between the spectrum of 
odd-parity Proca modes (with 5 = 0) and the spectrum 
of the massive scalar-field (s = 0) modes. In both cases, 
the exponent is 77 = Al + 5, although the constant of 
proportionality is larger for Proca-ficld perturbations. 

Figure [7] compares the ground state frequencies of the 
Proca field (I = 1, S = — 1) with those of the massive 
Dirac [2l| and scalar [H| fields. In the small-coupling 
limit, we see that, regardless of spin, all fields exhibit a 
hydrogenic spectrum, lj/h »1- (M/x) 2 /2. For low cou- 
plings, Mfj, < 0.4, the Proca field is more stable than the 
other fields, i.e., decays more slowly. At larger couplings, 
we see that the Proca field exhibits the largest binding 
energy (Rc(w//i— 1) « —0.14), although here the lifetime 
of the state is actually only a few black- hole light-crossing 
times. 
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Polarization for Bound States of Even-Parity Field 



Ground state energy levels for massive fields of spin = 0, 1/2 and 1 
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FIG. 5: Polarization state of even-parity bound states. The 
plot shows the complex number V, i.e. the ratio W(3)/ti(2) 
evaluated asymptotically [see Eq. (|46[) ]. as a function of the 
mass coupling M[i = . . . 0.8, for a selection of the lowest 
modes (I = 0, 1, 2) and overtones (n = 0, 1). The points show 
the values of V at Mfi = 0, 0.2, 0.4, 0.6 and 0.8. In the limit 
of vanishing mass, we find P — > — I for S = +1 modes, and 
V — !> / + 1 for S 1 = — 1 modes (where S is the spin projection). 
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FIG. 6: Comparison between the bound-state frequencies of 
the 1 = 1, n = modes for scalar (spin-0) and Proca (spin-1) 
fields. The imaginary part of the frequency exhibits the same 
dependence on the mass coupling M[i for both the scalar and 
longitudinal vector (odd-parity) bound states, being larger for 
the latter, while for the transverse vector (even-parity) states 
it grows with distinct powers of Mfi, being larger (smaller) 
for S = -l (S = +1) [see Eq. f5ty]. 
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FIG. 7: Ground state energy levels of massive fields of spin 
s = (scalar), s = 1/2 (Dirac), and s = 1 (Proca) on the 
Schwarzschild spacetime. The upper plot shows the real part 
of frequency Re(w//x) as a function of the mass coupling Mfj, 
and the lower plot shows (the absolute value of) the imaginary 
part of frequency Im(o;//i), which sets the decay rate, on a 
logarithmic scale. 



V. ANALYTICAL RESULTS FOR SMALL-MASS 
COUPLING 

To better understand the numerical results obtained in 
the previous section, we now study the equations for the 
odd- and even-parity modes in the limit of small mass 
(i and small frequency oj. The basic approach is to de- 
rive separate solutions in the near-horizon and far-field 
regions, defined by -C I and x 3> 1, respectively, 

where x = r/2M — 1, and then match them in their 
common domain of validity, 1 < i « As this re- 

quires lo I, we expect this to give a better approxima- 
tion to the numerical results for bound states, for which 
u ~ (j < 1, than for the QN modes, where u) > 0(0.1) 
even for (i 0. This approach will nevertheless give us 
a better insight into the analytical form of the solutions 
and the massless limit, being complementary to the more 
precise numerical analysis performed earlier. 

We start by writing the equations for the perturbations 
in terms of the dimcnsionless variable x, yielding for the 
odd-parity modes: 



[x 2 (x + ifdl + x(x + l)d x + V{x)]u (i) = 



(51) 
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where V(x) = 4u> 2 (x + l) 4 - A^i 2 x(x + l) 3 - X 2 x(x + 1), 
A 2 = 1(1 + 1) and both uj and [i are given in units of the 
inverse black hole mass M _1 . For the even-parity modes, 
one finds the coupled second-order differential equations: 

[x 2 (x + l) 2 d 2 x + x(x + l)d x + V(x)] u {2) = 

= x(2x - 1)(U( 2 ) - U( 3 )) , 

[x 2 (x + l) 2 d 2 ] +x(x + l)d x + V{x)] u (3) = 

= -2X 2 x(x + l)u {2) , (52) 

As we will see below, it will be convenient to write this 
system in terms of ip(x), defined in Eq. (|25j) . which up to 
a constant rescaling may be written as 



re) 



xU(3) 



X u 



(2) 



X + 1 



(53) 



Replacing U( 2 ) by this function, one obtains 
[x 2 (x + lfd 2 x + x(x + l)d x + V(x)] ij) = 

= 4:fi 2 x(x + l)u( 3) , 
[x 2 (x + lfdl+x{x + l)(2x + l)d x + V(x)]u {3) = 

= 2x(x + 1) 2 V • (54) 

This explicitly shows that the equation for ip decouples 
in the massless case, giving the physical 'vector' mode 
solutions described earlier, whereas for ip = one obtains 
a decoupled equation for U( 3 ) that yields the unphysical 
'scalar' modes. Let us now analyze the behaviour of the 
odd- and even-parity solutions in more detail. 

A. Odd-parity modes 



In the near-region, Eq. (fBTj) reduces to 
[x 2 (x + l) 2 d 2 x + x(x + l)d x + 4uj 2 - X 2 x(x + 1)] w (4) 

= , (55) 

which has a general solution given in terms of hypergeo- 
mctric functions 

u ™ ar = A {4) x- 2iu (x + l) 1+5 x 



2 F 1 (-l - 2iu + S, I + 1 - 2iw + S, 1 - Aiu, -x) + 

+±S(4)X [X + 1) X 

2 Fi (-1 + 2iuj + 5, 1 + 1 + 2iu + 6, 1 + Aiur, -x) , 

(56) 



where 8 = vT — 4oj 2 . It can be easily seen in terms of 
the tortoise coordinate r* that ingoing solutions at the 
horizon require setting = 0, and using the asymp- 
totic properties of the hypergeometric function 67[ we 
can derive the i> 1 form of the near-region solution 



vffi r ^ A (4) r[l - 4iw] x 



T[2l + 1] 



T[l + 1 - 2iuj + 5}T[l + l-2iuj-5] 

r[-2/-i] 

T[-l - 2iu + 5]T[-l - 2iuj - 6} 



.1+1 



(57) 



In the far-field region, Eq. (|51[) can be written as 

[x 2 d 2 x - 4q 2 x 2 + Aqvx - A 2 ]it (4) = , (58) 

where q and v were defined earlier [see Eq. (|30|) ]. and 
the general solution can be written in terms of confluent 
hypergeometric functions 

= e- z ^ 2 [C {4) z l+1 M(l + l-^2l + 2,z) + 

+ D {A) z~ l M(-l-v, -21, z)] , (59) 

where z = Aqx. For bound states, the linear combination 
which is regular at infinity corresponds to 

u bound = c {4)e - z/2 z l+1 U{l + 1 - v, 21 + 2, z) , (60) 
which for z <C 1 takes the form 16711 



bound 
X (4) 



(4) 



(4?) 



l + lr, 



J + l 



sin (21 + 2)tt 
(4 g )- 2 '- 1 
T[l + 1 - v]T[-2l\ 



T[-l - v]T[2l + 2} 



(61) 



Thus, the near- and far-region solutions yield the same 
power-law behaviour in the intermediate region, and one 
can equate the associated coefficients to get the matching 
condition 



T[-l - v]T[2l + 2} 
T[l + l- v]r[-2l] 

r[i + i- 



-(4q) 



21- 



-ir[-2Z-i] 



2iuj 



T[2l + 1] 
5} F[l + 1 - 2iu - 8} 



T[-l-2iuj + S] T[-l-2iu-8] 



(62) 



This condition can be solved in a similar way to [42j, 
taking into account that the left hand side vanishes to 
leading order for q <C 1, corresponding to the poles of 
T[l + 1 — v\. These are given by v = I + 1 + n for a 
non-negative integer n, yielding to lowest order in Mfj, a 
spectrum of Hydrogen-like bound states, as in the scalar- 
field case, and in agreement with our numerical results 



/' 



1 



2(i + i + n y 



(63) 



We may then expand both the left and right-hand sides of 
the matching condition about this value to get the next- 
to-leading order (NLO) correction to the spectrum. This 
is not as straightforward as for scalar-field perturbations, 
given that there is an uncancelled pole in one of the T 
functions. This can be overcome by taking to lowest order 
T[l + l-2iuj + S} ~ r[/+2] and F[-l-2iuj+S] ~ r[-Z + l], 
which is in fact consistent with the approximations used 
in the near-region. We may proceed as in [42j to obtain 
the next-to-leading order correction, valid for I > 1 



8uj 



4 2/ + V' +5 (2/ + 1 + n)\ (I + l)\(l - 1)! 



it 



■ (/ + l + n) 2 '+ 4 
/ 

x (1 + 2iw) Y[ (k 2 - 1 - Mw) , 

k=l 



[(2/)!(2/ + l)!] 2 



(64) 
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where u> takes the leading order value in Eq. (|63|). From 
this one can extract the imaginary part of the bound- 
state frequency, which to leading order grows as fi 4l+6 
as for scalar-field perturbations, although with a differ- 
ent coefficient. For example, for the lowest-lying dipole 
mode, I = l,n = 0, we obtain ujj ~ /i 10 /3, which is 
twice the value obtained for scalar-field perturbations 
[4(|[42|. In this sense, we classify the odd-parity modes as 
longitudinal bound states, which behave like scalar-field 
perturbations far from the black hole but have a vec- 
torlike near-horizon behaviour, which makes the lowest- 
lying modes decay more rapidly. We illustrate this be- 
haviour in Fig. [8j where one can see that the numerical 
curve approaches the matching result for small masses. 



Decay rate of odd-parity 1=1 bound state 
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FIG. 8: Comparison between the numerical and analytical 
results for the imaginary part of the bound-state frequency 
for the odd-parity I = 1, n = mode as a function of the 
mass coupling M[i. The solid [red] line shows numerical data, 
and the dashed [blue] line shows the analytical approximation 
Im(w//i) « -\{Mp,f. 



neglect v and obtain the QN mode condition 
1 = (4c?) 



2i+1 r[-2?] r[z + i] r[-2Z-i] 



r[-i] r[2/ + 2] r[2/ + i] 

T[l + 1- 2iuj + 5]T[l + l- 2iu - 5} 
T[-l - 2iu + 8} r[-/ - 2iuj - 5] 



(66) 



We may then proceed as for bound states to approxi- 
mately eliminate the uncancelled V pole, obtaining to 
leading order, for I > 1 



1 = -i(4Vw 2 -M 2 ) 2i+1 (I + 1)1(1 - 1)1 

2 I 

(l + 2icj) Y[{k 2 - 1 -Uu) 
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(2i)!(2i+l)l 



k=l 



(67) 



The solutions of this equation in the lower half of the 
complex plane then yield the spectrum of QN 

modes for odd-parity perturbations. For example, for 
1 = 1, solving this to leading order in u> and /i gives 



odd 
(1,0) 



0.515(1 + 0.364(1 + i) f i 2 



(68) 



This result overestimates the real and imaginary part ob- 
tained numerically, which was expected given that the 
matching procedure only holds for wCi, but it neverthe- 
less gives the correct qualitative /i-dcpendcncc of the QN 
mode frequencies. Notice also that this does not agree 
with the matching analysis for masslcss vector fields per- 
formed in [68} using the Teukolsky equation rather than 
the spin-1 Regge- Wheeler equation. However, this anal- 
ysis also departs significantly from the numerical data, 
which again is symptomatic of the failure of the approx- 
imations involved. 



B. Even- parity modes 



For the QN modes, the coefficients Cm and Dm in 
Eq. (|59|) can be obtained by taking into account that 
M(a, b, z) ~ 1 for small z and then matching them with 
the corresponding coefficients in Eq. ([57)) . Asymptoti- 
cally, the far-region solution takes the form 



QN 

V) 



J / 2 (-l)-' y+i+ V X 



a 



T[2Z + 2] 



rp + i + v] 



D 



(-1) 



r[-2*] 

T\-l + u] 



e z ' 2 z- v 



C, 



(4) 



T[2l + 2] 

r[i + i- v] 



D 



(4) 



r[-2/] 



(65) 



Recalling that z = ^Iqx 



\x , it is easy to see 



that the first term corresponds to ingoing waves from 
infinity and should be set to zero. One should notice that 
for QN modes one expects y<!in the limit /i«w<I. 
as opposed to the case of bound states where v takes 
integer values as seen above. We may thus consistently 



The problem becomes more involved for the even- 
parity modes, where one has a system of coupled dif- 
ferential equations and analytical results are harder to 
extract using the matching procedure described above. 
In the far-region, however, the system may be diagonal- 
ized and a general solution can be written as 



2,3) (z) = C (2,3) U (S)( Z ) , 



(69) 



s=±i 



where urs) satisfies the confluent hypergeometric equa- 
tion: 

[x 2 d 2 x - Aq 2 x 2 + Aqisx - j(j + l)}u (s) - , (70) 

with j = I + S as defined earlier, and cCL/ci^ = —I, 
c (3) / c (2) = ^ + 1 • The general solution is analogous to the 
odd-parity far-region solution, yielding: 

U (S) = e-^ 2 [C (S) z j+1 M(j + l-v,2j + 2,z) + 

+ D {S) z-iM(-j-v,-2j,z)] , (71) 
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with bound states corresponding to the confluent hyper- 
geometric function U (j + 1 — v, 2j + 2, z). 

In the near region, one could hope that the solutions 
would converge to the massless case expressions, given 
that the mass term can be neglected in this limit. How- 
ever, writing the system of equations as in Eq. (|54j) . one 
explicitly sees that the effects of the vector field mass may 
only be neglected if the functions ip and are of com- 
parable magnitude, which a priori is not necessarily the 
case. We may nevertheless investigate the form of these 
solutions neglecting the mass term and check whether a 
matching procedure can be used. 

The generic solution in the near-region for ip is then 
identical to that obtained for U(4), yielding for ingoing 
boundary conditions at the horizon: 



iV* (z + 1) 



1+5 



iF x {-l - 2iu + 6, 1 + 1 - 2iu> + 6, 1 - -x) , 



constituting the 'vector' solution described earlier. One 
may also derive the form of the 'scalar' solutions in the 
near-region by setting ip = and solving the decoupled 
equation for tt(3), which may also be written in terms of 
a hypergeometric function: 



'(3) 



(x + 1) 5 



x 2 Fi(-l,l + l,l-4ioj, 



(73) 



As expected, this corresponds to the solution for scalar- 
field perturbations in the near- region. For large x, these 
solutions behave like 1671: 



A^l - Aico] x 

r[2i + i] 



l (3) 



r[Z + 1 - 2iu + 6}T[l + 1-2ilo-6] 

r[-2Z-i] _ ; 

r[-Z - 2iuj + 5]T[-l - 2iu> -6} X 
A(3)T[1 - Mw] x 

r[2Z + i] , , 



T[l + 1}T[1 + 1 - 4iw] 
r[-2/-l] 
T[-l]T[-l - 4ioj] 



(74) 



One may then use the far-region solutions in Eq. (|6T))) 
to compute these quantities in the far-region and take the 
limit z <C 1 as for the odd-parity modes. It is then easy 
to see that matching is only possible for v <C I which, 
from our analysis of the odd-parity modes, is expected 
to be the case for QN modes but not for bound states, 
where v takes positive integer values. 

That a matching between the (massless) near-region 
solutions and those in the far-field region is not possi- 
ble for bound states is not completely unexpected, since 
these modes have no electromagnetic analogue as mass- 
less waves cannot be bound in a gravitational field. Fur- 
thermore, the form of the far-region solutions actually 



suggests a nontrivial mixing between the 'vector' and 
'scalar' solutions of the massless case that cannot be de- 
termined using this approach. 

Nevertheless, were the near-region solutions fully 
known, one could construct a matching condition analo- 
gous to Eq. (|62p for each of the far-region solutions U(±), 
labeled by j = I ± 1. Thus, one expects the bound- 
state spectrum to be given, to lowest order, by the poles 
of T[j + 1 — v], which yield the hydrogenic spectrum in 
Eq. (|49|) . in agreement with the numerical results. Recall 
that in the limit Mfj, <C 1 we had obtained V = — I, l + l, 
which is consistent with pure U(±) bound-state solutions 
far from the black-hole horizon. The polarization data 
suggests, however, that as Mfj, increases the mixing be- 
tween these solutions becomes larger. 

The imaginary part of the even-parity bound-state 
modes depends, unfortunately, on the particular form of 
the near-region solution. One would naively expect it to 



(72) S row likc M 



4j+6 



41+4S+6 



in analogy with the odd- 



parity case, but as discussed in the previous section one 
can infer an additional factor of n~ 2S from the numerical 
data, which is clearly suggestive of a nontrivial mixing of 
the 'vector' and 'scalar' near-region solutions. 

On the other hand, for QN modes one may in principle 
take the limit v <C I, where the confluent hypergeometric 
functions in Eq. I 
Bessel functions 



can be written in terms of modified 



M(n + l,2n + 2, z) = 

r[n + l/2] e '/'(f 



-n-l/2 



W/a(|) (75) 



and we may write the general solution as: 



far r~ 

U (S) = v z 



C {s) I a+s {z/2) + b [S) I^ s {z/2) 



(76) 



where a = 1 + 1/2. One may then use the relations 
between adjacent Bessel functions (67} to show that 







far 



C v I a (z/2) + D v I- a (z/2) = 
C v z l+1 M(l + l 1 2l + 2,z) + 



+ D v z~ l M(-l,-2l,z) 



(77) 



where, setting = 1 without loss of generality, we have 



C v = -IC W + + 1)(7 ( _) , 

Dy = -lb { + ) + + 1)%) . 



(78) 



The corresponding constants Cv and Dy can be obtained 
from these using Eq. (|75|) . This shows that, within this 
approximation, ip has the same form as the odd-parity 
function it( 4 ) both in the near-horizon and asymptotically 
flat regions and should hence yield the same QN mode 
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spectrum. As our numerical analysis shows, this holds 
only in the massless limit, so that this approach fails 
to describe the broken degeneracy between the even- and 
odd-parity 'vector' QN modes for finite fi. It nevertheless 
illustrates how the two j = Z±l solutions combine to form 
'vector' states in the far-field region for small v. 

One may also obtain the 'scalar' QN modes corre- 
sponding to solutions with ip = 0. In the far region 
this implies from Eqs. ([77]) and (78]) that C , (+) /C'(_)- = 
Dr + \/Dr_\ = (I + and we may write 



far 
l (3) ~ 



2(l + l)(2l + l) Z - 1/2 
C ( _ ) I m/2 ( 2 /2) + J D ( _ ) /_ i _ 1/2 ( Z /2) 



:/2 



C\ 3) z l M(l + l,2l + 2,z) + 



+ D {3) z- l - 1 M(-l~2l,z) 



(79) 



We may then use the asymptotic properties of the mod- 
ified Bessel functions [63] in this case to show that 
u {3)/ u (2) ~* 0, in agreement with the numerical polar- 
ization data in the massless limit. Taking the inverse 
ratios C'( + )/C , (_) = D( + )/Dr_) = 1/(1 + 1), we obtain 
u (3)/ u (2) 1 asymptotically, which correspond to the 
'vector' even-parity QN modes as obtained numerically. 

We may derive a matching condition for the 'scalar' QN 
modes in a similar fashion to the odd-parity case, by first 
matching the constants C( 3 ) and -D(3) to the correspond- 
ing coefficients in Eq. (1741) and then imposing outgoing 
waves at infinity. From Eq. (I79p . we obtain asymptoti- 
cally 



far 
L (3) 



z/2 



"(-I) 



l+l 



(3) 
,z/2 



T[2l + 2} 



)-^(3) 



r[-2*] 
r[-i] 



(3)" 



T[2l + 2} 



D 



T\-2l] 



(3) 



n-i] 



(80) 



Setting the coefficient of the first term to zero, we ob- 
tain the 'scalar' matching condition 



1 = (4g) 

r 



a+iEhgj] r[ 



i] r[-2i-i]v[i + i] 



v[-i] r[2Z + 2] r[2/ + i] r[ 

1 - 4du] 



-1} 



T[-l - Aiuj] 
This can be further simplified to yield 



(81) 



(82) 

For comparison with the 'vector' QN modes, we solved 
this equation to leading order for the first excited state, 
yielding: 

even^S 



'(1,0) 



0.612(1 + 0.306(1 + , (83) 



which as before overestimates the numerical result but 
yields the correct qualitative [i dependence, also showing 
that the 'scalar' QN modes have larger frequencies than 
the corresponding 'vector' states, in agreement with our 
earlier results. 

Hence, although the analytical matching procedure 
cannot really replace the numerical analysis in terms of 
quantitative results, it illustrates the rich and nontriv- 
ial structure of the massive vector field perturbations, in 
particular that of the even-parity modes. This interesting 
structure is a consequence of the different spin structure 
near and far from the black-hole horizon, leading to a 
nontrivial interplay between spin and orbital angular mo- 
mentum. In particular, whereas for QN modes one finds 
'scalar' and 'vector' states at infinity as in the massless 
case, these two components are nontrivially mixed for 
bound states, which have no massless counterpart. 



VI. CONCLUSION AND FUTURE PROSPECTS 

As we have found, massive vector fields exhibit 
an extremely rich spectrum of perturbations on the 
Schwarzschild spacetime, due to both the nonvanishing 
mass and the spin-1 nature of the field. Our results show 
that the 'vector' and 'scalar' solutions describing elec- 
tromagnetic perturbations on this geometry, the latter 
being unphysical gauge degrees of freedom in electromag- 
netism, mix in a nontrivial fashion in the presence of a 
nonzero mass, even for small-mass coupling Mfj,. 

For the electromagnetic field (/i = 0), (i) the even- and 
odd-parity 'vector' modes are governed by the same dy- 
namical equation, and hence their QN frequency spectra 
are degenerate, and (ii) the 'scalar' degree of freedom 
corresponds to a 'pure-gauge' mode. By contrast, in the 
Proca case (/j ^ 0), as a consequence of the breaking of 
the electric-magnetic duality, we find that (i) the odd- 
and even-parity 'vector' modes become dynamically dis- 
tinct, and (ii) the even-parity 'scalar' mode acquires a 
physical significance, since the Proca field has no gauge 
freedom (A? = 0). We showed in Sec. Ill Bl that the odd- 
parity part of the field is governed by a single equation, 
whereas the even-parity part of the field is determined 
by a coupled pair of equations. 

The Proca field exhibits quasibound states, i.e. solu- 
tions which can be localized within the vicinity of the 
black-hole horizon and which are absent in the Maxwell 
case. Quasibound states on the Schwarzschild space- 
time have complex frequencies, with an imaginary part 
corresponding to the decay rate (as flux is absorbed 
by the horizon). Our numerical and analytical studies 
of the bound-state spectra reveal an interplay between 
the 'vector' and 'scalar' solutions, which is reminiscent 
of a spin-orbit coupling between the field's proper spin 
and the orbital angular momentum of each multipole. 
States may be labeled by their total angular momentum 
j = I + S, as measured by an asymptotic observer. We 
thus find 'electric' (even-parity) transverse states, with 
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j = Z±l, and 'magnetic' (odd-parity) longitudinal states, 
with j = I, in both cases yielding a hydrogenic spectrum 
for Mfi <C 1 labeled by the 'principal quantum number' 
N = j + 1 + n for non-negative integers n. While this 
agrees with earlier studies for both the monopole [3l[ 
and large multipoles [l7j |. we find decay times which are 
parametrically different for each type of mode, as op- 
posed to [l7j where a common behaviour for small-mass 
coupling Im(w//i) oc (M fi) 4l+5 was found for all spin-j 
states, albeit with different coefficients. Although some 
further analytical insight is required to better understand 
this behaviour, from our numerical results we can infer 
a power-law behaviour Im(o;//x) oc (M /j,) 4:l+2S+5 in the 
same limit. 

The fact that the bound-state decay rate is sensitive to 
its spin in addition to its orbital angular momentum can 
be understood in simple physical terms in the 'antitunnel- 
ing' picture devised in |48[. Bound states are localized in 
a potential well whose depth depends on the field's mass 
and total spin and which is separated from the black- 
hole horizon by a finite angular momentum barrier. The 
rate at which the black hole attenuates the mode's wave- 
function is then determined by the height of this barrier, 
which is controlled by the total angular momentum of the 
state. In particular, the states with the smallest angu- 
lar momentum j — I — 1 for each multipole are absorbed 
more quickly by the black hole and hence exhibit a faster 
decay rate. For higher multipoles, we have nevertheless 
that j ~ I ^> 1, which justifies the results obtained in 
13. 

These results also suggest that a similar behaviour 
should be observed for rotating black holes, where wave 
modes with u> < mQ are amplified rather than damped 
by superradiant scattering in the Kerr ergoregion, with 
—l<m<l denoting the azimuthal angular momentum 
projection and f2 the black hole's rotational frequency. 
When such states are bound to the black hole, multi- 
ple wave scatterings will amplify the corresponding wave 
function and consequently exponentially increase the as- 
sociated particle number, giving rise to the so-called black 
hole bomb effect [3(| . The total angular momentum bar- 
rier should then also determine the overlap of each mode 
with the ergoregion and we thus expect the j = I — 1 
bound states of the Proca field to exhibit a parametri- 
cally faster instability rate in this case as well. 

A rigorous analysis of the massive vector field bound 
states on the Kerr spacetime poses, however, an ex- 
tremely challenging problem even from the numerical 
point of view as the field equations do not seem to admit 
separable solutions in this geometry. In addition, parity 
invariance is broken by the black hole's rotation, which 
may give rise to nontrivial mixings between the odd- and 
even-parity states, namely in the extremal case, where 
we expect the superradiant instability to be strongest 
13, Hi| . One might expect these effects to become sub- 
dominant for slowly-rotating black holes, and a prelimi- 
nary analysis in this limit seems to confirm our physical 
intuition, with the lowest spin states exhibiting a para- 



metrically larger instability, although it remains unclear 
whether it is consistent to study superradiant modes in 
this case. 

Although a comprehensive study of Proca perturba- 
tions on the Kerr spacetime is beyond the scope of this 
work, the physical picture derived from our results for 
nonrotating geometries suggests that Kerr superradiance 
may be relevant for probing the existence of ultralight 
hidden U(l) vector fields in string compactifications, as 
described in the introduction. In particular, one expects 
this effect to be more pronounced for hidden photons 
than for axionlike fields, given the existence of states 
with a lower angular momentum barrier due to the above 
mentioned spin-orbit coupling. In addition, the lowest- 
lying odd-parity longitudinal bound states of the Proca 
field exhibit a faster decay rate than the corresponding 
scalar-field modes in the Schwarzschild case as, despite 
yielding similar angular momentum barriers, higher-spin 
waves decay more rapidly in this case, a well-known re- 
sult for massless fields |38j . 

This suggests that, in the case where axionlike and hid- 
den photons of similar masses coexist, superradiant scat- 
tering will amplify the latter bound states more quickly. 
This will in turn inhibit the formation of axionlike bound 
states, given that the hidden photon cloud breaks the 
rotational symmetry of the system and thus suppresses 
multiple scatterings in states with distinct quantum num- 
bers, as described in [49[. The cloud may later be de- 
pleted by self-interactions and the progressive shutdown 
of the mode instability due to the decreasing black hole 
mass, so that other states may be amplified, possibly 
leading to an interplay between spin-1 and spin-0 states 
with interesting phenomenological consequences. Fur- 
thermore, similarly to the axion case, massive hidden vec- 
tor field clouds around astrophysical black holes should 
also lead to phenomena such as gravitational waves and 
possibly 'bosenova-like' emission, while the generic mix- 
ing between hidden and visible photons could yield ex- 
citing novel signatures j48l - l50j . 

The study of massive vector field perturbations in 
black hole spacetimes is thus an important problem from 
both the black hole stability and phenomenological per- 
spectives. We hope in the future to further develop and 
extend our analytical and numerical methods to better 
understand the behaviour of massive higher-spin pertur- 
bations on the Schwarzschild spacetime, as well as for 
more generic black hole geometries with angular momen- 
tum and charge. 
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